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SPACE-TIME  ANALYSIS  OF  THE  FLOW-INDUCED 
VIBRATION  OF  A  SOLID  CONTAINING  AN 
EMBEDDED  CYLINDRICAL  SHELL 

1.  INTRODUCTION 

The  analysis  of  solid  bodies  with  embedded  cylindrical  shells  that  contain  transported 
fluids  is  not  only  pertinent  to  such  mechanical  systems  as  buried  steam  pipes,  sewers,  and  heat 
exchangers,  but  also  to  biological  systems,  such  as  the  arteries  of  the  human  body.  The 
response  of  the  solid  body  is  dependent  on  the  forces  that  the  transported  fluid  exert  on  the 
surrounding  medium.  These  forces  typically  cause  wave  motion  in  the  medium,  which,  in 
turn,  creates  a  displacement.  Understanding  the  response  of  the  medium  is  particularly 
important  when  the  displacement  and  wave  motion  have  undesirable  effects.  The  response  of  a 
free  surface  of  the  medium  may  also  be  used  to  detect  and  localize  the  source  (energy) 
emissions.  The  work  documented  in  this  report  develops  the  theory  and  provides  the 
experimental  data  to  support  detection  and  localization  of  flow  in  an  unblocked  shell 
surrounded  by  a  lossy  medium.  Although  this  type  of  detection  may  be  useful  to  surveyors 
trying  to  locate  buried  pipes,  its  primary  intent  is  to  locate  partially  blocked  coronary  arteries  in 
humans. 

In  section  2,  a  model  of  a  solid  with  an  embedded  tube  containing  a  moving  fluid  is 
derived  based  on  the  wavenumber  (space)  and  frequency  (time)  content  of  the  wave  motion. 
The  model  is  evaluated  at  the  siuface  of  the  solid,  which  corresponds  to  the  position  of  a 
beamformed,  linear  array.  The  model  is  developed  using  the  equations  of  elasticity  in 
cylindrical  coordinates  coupled  to  a  line  source  that  represents  the  fluid.  The  medium  is 
modeled  as  a  semi-infinite  solid  that  admits  propagating  shear  and  compressional  wave  energy. 
The  filtering  effects  of  the  array  and  sensors  are  also  included  in  the  model. 

In  section  3,  the  closed-form  solution  of  the  model  is  compared  to  experimental  data  for 
validation.  The  experiment  consists  of  a  solid  Hexcel  195-RE  urethane  block  that  contains  an 
embedded  surgical  tube  (figure  1).  This  block  is  an  isotropic  solid  free  of  bone,  tissue,  and 
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other  biological  substances  whose  behavior  cannot  be  easily  resolved.  The  tube  has  an  inner 
radius  of  0.00159  m  (1/16  in.),  which  is  the  approximate  size  of  a  human  coronary  artery.  The 
distance  from  the  surface  of  the  solid  to  the  tube  is  0.0191  m  (3/4  in.),  which  is  approximately 
the  minimum  distance  from  a  location  on  the  third  intercostal  space  to  the  nearest  coronary 
artery  in  a  human  chest  (Norton,  1995).  Water  is  gravity  fed  into  the  surgical  tube  from  a 
holding  tank  and  captured  below  the  block  in  a  calibrated  flask  (figure  2).  This  technique 
permits  measurement  of  the  flow  velocity.  The  moving  fluid  transmits  a  relatively  small 
normal  force  into  the  embedded  tube,  which,  in  turn,  transmits  energy  into  the  solid.  This 
energy  propagates  to  the  surface  of  the  block,  primarily  in  the  form  of  shear  waves  where  it 
produces  a  local  displacement. 

The  data  in  this  experiment  were  collected  with  an  eight-element  array  (figure  3) 
designed  to  detect  the  extremely  small  displacements  that  are  typical  of  experiments  with  a  low- 
level  energy  source.  Each  element  of  the  array  consists  of  a  piece  of  stiff  plastic  mounted  on  a 
strip  of  polyvinylidene  fluoride  (PVDF).  The  PVDF  is  a  piezoelectric  polymer  that  acts  as  a 
strain  gage  when  the  plastic  is  displaced.  Array  technology  improves  the  signal-to-noise  ratio 
over  single-point  measurements  and  provides  directionality  to  the  detection  process. 

In  section  4,  the  dynamic  behavior  of  the  model  is  discussed,  and  the  differences 
between  the  model  and  the  experimental  data  are  examined. 

The  model  and  experiment  developed  here  allow  a  controlled  laboratory  setting  in 
which  array  technology  is  used  to  conduct  tests  of  arterial  flow.  This  effort  provides  a 
foundation  for  future  investigations  that  would  address  the  vibration  of  shells  that  include 
partial  blockages,  nonisotropic  blocks,  and  focused  optimal  processing  (beamforming) 
detection  techniques. 
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Calibrated 
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Figure  2.  Laboratory  Configuration 
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Figure  3.  The  Eight-Element  Array 
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S(k„Q)o)=  \G(kf  \r{k,(OQf  \Ptbi{k,(OQf  dk  ,  (1) 

— <o 

A 

where  ^(k^,0}Q)  is  the  beamformed  displacement  wavenumber- frequency  spectrum 
corresponding  to  the  measured  output  at  steered  wavenumber  k^  (rad/m)  and  fixed  frequency 
C0()  (rad/s),  H(k-kg)  is  the  steered  response  pattern  of  the  array,  G(k)  is  the  response  pattern 
of  each  individual  sensor  in  the  array,  T(k,C0Q)  is  the  transfer  function  between  the 
displacement  of  the  solid  at  the  surface  and  the  normal  pressure  of  the  turbulent  boundary  layer 
of  the  fluid  acting  on  the  embedded  cylindrical  shell,  Ptbi(f^,0)Q)  is  the  streamwise  turbulent 
boundary  layer  pressme  in  the  fluid,  and  k  is  the  wavenumber  (rad/m).  Each  of  these  terms  is 
developed  below. 


2.1.  THE  STEERED  ARRAY  RESPONSE 

The  normalized  wavenumber  response  of  a  linear  array  of  N  discrete  sensors  can  be 
modeled  as 


N 


.iik-ks)Xfi 


N 

«=1 


(2a) 


where  is  the  position  of  the  center  of  sensor  n,  and  is  a  spatial  weighting  function  (e.g., 
Taylor  or  Hanning).  For  this  analysis,  is  unity  (Boxcar  window).  For  a  linearly  spaced 
array,  the  wavenumber  response  with  uniform  shading  is 


H{k-k,)  = 


(2b) 
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where  d  is  the  center-to-center  separation  distance  of  adjacent  sensors  (m).  Equations  (2a)  and 
(2b)  assume  that  the  sensors  are  responding  as  point  sensors.  Figure  4  is  a  plot  of  equation 
(2b)  versus  wavenumber.  The  sensor  separation  distance  {d)  is  0.0127  m  (0.5  in.)  and  the 
number  of  sensors  in  the  array  (AO  is  eight.  The  main  beam  of  the  array  has  been  steered  {k^ 
to  150  rad/m.  The  dashed  line  corresponds  to  the  unaliased  region  of  the  array,  which  is  equal 
to  ±K Id  or  ± 247.4  rad/m. 


Unaliased 


Figure  4.  The  Theoretical  Beam  Response  of  an  Eight-Element 
Array  Steered  to  a  Wavenumber  of  150  rad/m 


7 


TR  11,073 


2.2.  THE  INDIVIDUAL  SENSOR  RESPONSE 

The  response  of  each  individual  sensor  is  modeled  as 


(3) 


where  L  is  the  length  of  the  sensor  (m).  Figure  5  is  a  plot  of  equation  (3)  versus  wavenumber. 
The  sensor  length  (L)  is  0.0127  m  (0.5  in.).  The  dashed  line  corresponds  to  the  unaliased 
wavenumber  region  of  the  array  described  in  section  2.1. 


Unaliased 
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2.3.  THE  TRANSFER  FUNCTION 

The  transfer  function  between  the  vertical  displacement  of  the  solid  at  the  top  surface 
and  the  turbulent  boundary  layer  pressure  of  the  internal  fluid  is  derived  by  assuming  that  the 
system  is  composed  of  two  concentric  cylindrical,  linear,  isotropic  mediums,  each  governed  by 
the  equation  (Timoshenko  and  Goodier,  1934) 


(4) 


where  p  is  the  density;  X  and  p  are  the  Lam6  constants;  t  is  time;  •  denotes  a  vector  dot 
product;  u  is  the  cylindrical  coordinate  displacement  vector  expressed  as 


u  =  < 


Urir,e,z,t) 

UQir,d,z,t) 

U2{r,d,z,t) 


(5) 


with  subscript  r  denoting  the  radial  direction,  6  denoting  the  angular  direction,  and  z  denoting 
the  axial  direction;  V  is  the  gradient  vector  differential  operator  written  in  cylindrical 
coordinates  as  (Potter,  1978) 


yj  d  .  1  .  d  . 


(6) 


with  ij.  denoting  the  unit  vector  in  the  r-direction,  denoting  the  unit  vector  in  the  6-direction, 
and  denoting  the  unit  vector  in  the  z-direction;  V  is  the  three-dimensional  Laplace  operator 


operating  on  vector  u  as 


2„  Mr  2  dug] 


v^u=  vx--^- 


r2  r^der 


\lr  + 


+  +  (7) 


r\ 

with  V  operating  on  scalar  u  as 


-  V .  Vu,  -  '  ^  ]  ,  1  ^\.ej  ,  . 


and  the  term  V  •  u  is  called  the  divergence  and  is  equal  to 


(8) 


„  dUr  1  dua  du,  Ur 
V*u  =  -X-t- — X  + 

dr  r  36  dz  r 


(9) 


9 


TR  11,073 


The  pressure  effects  of  the  fluid  acting  on  the  interior  surface  of  the  inner  cylinder  are  modeled 
as  structural  loads  later  in  this  section.  The  modeled  geometry  and  the  coordinate  system  of  the 
cylinders  are  shown  in  figure  6.  Note  that  p,  X,  and  p  are  the  properties  of  each  specific 
material.  When  the  modeled  region  is  the  embedded  cylinder  {a<r<b),  the  material 
properties  of  the  surgical  tube  are  used.  When  the  modeled  region  is  the  solid  block 
(ft  <  r  <  c) ,  the  material  properties  of  the  urethane  are  used.  The  flat  surface  of  the  block  is 
modeled  as  a  cylinder,  which  results  in  this  surface  being  approximated  with  a  curve.  The 
effects  of  this  approximation,  which  is  necessary  to  obtain  a  closed  form  solution  of  the 
displacement  field  on  the  block,  are  discussed  later. 

The  displacement  vector  u  is  written  as 

u  =  V^  + VxH  ,  (10) 

where  ^  is  a  dilatational  scalar  potential,  x  denotes  a  vector  cross  product,  and  H  is  an 
equivoluminal  vector  potential  expressed  as 


H  = 


Hr(r,d,z,t) 

^  HQ{r,e,zj)  > 

H  (r,  9,  z,t) 


(11) 


Expanding  equation  (10)  and  breaking  the  displacement  vector  into  its  individual  terms  yields 


#  1  atf,  3H„ 

'  dr  r  de  A  ’ 

I  dill  dH,  dH, 
and 


(12) 

(13) 


Mz 


dz  r  dr 


IdHr 
r  de  ' 


(14) 
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Figure  6.  Modeled  Geometry  With  Coordinate  System 
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Equation  (10)  is  next  inserted  into  equation  (4),  which  results  in 


and 


d^(t> 

dP- 


(15) 


cJV^H  = 


d^H 

dP 


(16) 


on  both  the  inner  and  the  outer  cylinders.  The  constants  and  are  the  complex  dilatational 
and  shear  wave  speeds,  respectively,  and  are  determined  by 


and 


(17) 


— 


= 

p 


(18) 


The  relationship  of  the  Lame  constants  to  the  compressional  and  shear  moduli  is  shown  as 
Em 


(l  +  u)(l-2u) 


(19) 


and 


/i  =  G  = 


2(1 +  u)  ’ 


(20) 


where  E  is  the  complex  compressional  modulus  (N/m^),  G  is  the  complex  shear  modulus 
(N/m^),  and  v  is  the  Poisson's  ratio  of  the  material  (dimensionless). 

The  conditions  of  infinite  length,  axisymmetric  response  (n  =  0),  and  steady-state 
response  are  now  imposed,  allowing  the  scalar  and  vector  potential  to  be  written  as 


=  g(r)cos(nd)e^^e^‘^  =  g{r)P^e^^  , 

(21) 

Hr  =  hr  (r) sin(ne)e^^P^  =  0  , 

(22) 

HQ=hQ{r)cos{nd)P^P^  , 

(23) 
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and 


=  /i2(r)sin(rtd)e‘^e‘‘“  =  0  ,  (24) 

where  k  is  the  wavenumber  of  excitation,  co  is  the  frequency  of  excitation,  and  i  is  the  square 
root  of  -1.  Note  that  equations  (21)-(24)  are  valid  on  both  the  inner  and  the  outer  cylinders. 
For  axisymmetric  response,  the  equations  of  motion  are  dependent  only  on  the  scalar  potential 
<p  and  angular  contribution  Hg  of  the  vector  potential.  Additionally,  because  Hf.  =  0,H^  =  0, 
and  Hq  and  0  are  not  functions  of  6,  equation  (13)  becomes 

■<e  =  ^  =  o,  (25) 

where  ( )  denotes  any  function.  Inserting  equations  (21)-(24)  into  equations  (15)  and  (16) 
yields  the  following  four  distinct  wave  equations  on  both  the  inner  cylinder  and  the  outer 
cylinder: 


d  g(r)  ^  1  dg(r)  ^ 


dr^ 


r  dr 


j 


gir)  =  0  , 


d%{r)  ^  1  dh^jr)  ^ 


dr^ 


r  dr 


CO 


2  2 
r 


\Kir)  =  Q  , 


d\{r)  .  1  dhg{r) 


dr^ 


r  dr 


he(r)^0  , 


and 


(26) 

(27) 

(28) 


d\(r)  .  1  dh,ir)  , 
dr^  r  dr 


^  (0^ 


\ 


Mr)  =  0  . 


(29) 


The  solutions  to  equations  (26)  and  (28)  are  now  found  on  both  cylindrical  domains 
with  Bessel  functions.  No  solution  is  found  to  equations  (27)  and  (29)  because  they  do  not 
contribute  to  the  axisymmetric  response.  The  solution  to  equation  (26)  is 
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C\JQ{oc\r)  +  C2YQ{(Xir)  a<r<b 

g(r)  =  ^  ,  (30) 

CsJq  {0C2r)  +  Cglo  («2'‘)  b<r<c 

where  Jq  is  a  complex,  zero-order,  first-kind  standard  Bessel  function;  Yq  is  a  complex,  zero- 
order,  second-kind  standard  Bessel  function;  C^,  C2,  C5,  and  Cg  are  complex  constants 
(determined  below);  is  defined  on  the  inner  cylinder  as 


where  the  subscript  d\  refers  to  the  complex  dilatational  wave  speed  in  the  inner  cylinder;  and 
02  is  defined  on  the  outer  cylinder  as 


where  the  subscript  d2  refers  to  the  complex  dilatational  wave  speed  in  the  outer  cylinder.  The 
solution  to  equation  (28)  is 

C^MM  +  C^YiiM  a<r<b 

heir)  =  <  ,  (33) 

CiJiip2r)  +  CsYiHhr)  b<r<c 

where  /j  is  a  complex,  first-order,  first-kind  standard  Bessel  function;  Y^  is  a  complex,  first- 
order,  second-kind  standard  Bessel  function;  C3,  C4,  C7,  and  Cg  are  complex  constants 
(determined  below);  is  defined  on  the  inner  cylinder  as 


where  the  subscript  si  refers  to  the  complex  shear  wave  speed  in  the  inner  cylinder;  and  P2  is 
defined  on  the  outer  cylinder  as 
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where  the  subscript  ^2  refers  to  the  complex  shear  wave  speed  in  the  outer  cylinder. 

The  strains  are  related  to  the  displacements  in  an  axisymmetric  solid  by  (Timoshenko 
and  Goodier,  1934) 


and 


_ 


dr  ’ 


^ee  = 


_“r 


_ 


■'ZZ 


dz  ' 


(36) 

(37) 

(38) 


1 

£„  =- 
2 


(du,dUr 


(39) 


dr  dz 

The  relationship  between  the  displacements  (and  the  derivatives  of  the  displacements)  and  the 
potential  functions  g  and  hg  is  found  by  combining  equations  (12),  (14),  (21),  and  (23)  to 
produce 


and 


Ur  =( -  iichg  (r) , 


V  dr 


dUr  _ 


^  dr'^  dr  ) 


OT  V.  r  dr  ) 


dr 


.^dg{r)  ^\dhg{r)  heir)  ^d^heir) 


dr  r  dr 


dr^ 


Jkzeicot 


^  =  r  +  k\  (r)  V' 

dz  \  dr  J 


ikzgiox 


(40) 

(41) 

(42) 

(43) 

(44) 


(45) 
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The  displacements  and  external  forces  are  now  equated  by  use  of  the  stress-strain 
constitutive  equations  on  the  free  surfaces  of  the  cylinders  and  at  the  inner  and  outer  cylinder 
interface.  The  normal  stress,  strain,  and  radial  forces  acting  on  the  inner  cylinder  at  radius  a 
are  related  by 

(T pF  i^ci,  z,  t)  ( Aj  +  2^1  (fl)  z,  t)  +  ^^£qq  (fl,  z,  r)  +  AjCjz  —  Pi  »  (4ti) 

where  Orr(a,2,t)  is  the  normal  radial  stress,  err(ci,z,t)  is  the  normal  radial  strain,  eQQ(a,z,t) 
is  the  normal  circumferential  strain,  e^zia^z^t)  is  the  normal  longitudinal  strain,  Pi{a,i,t)  is  the 
internal  pressure  acting  on  the  cylinder  in  the  radial  direction,  and  the  subscript  1  denotes 
material  properties  of  the  inner  cylinder.  The  shear  stress,  strain,  and  longitudinal  forces  acting 
on  the  inner  cylinder  at  radius  a  are  related  by 

0  >  (47) 

where  Grz(.a,z,t)  is  the  shear  stress,  erz{a,z,t)  is  the  shear  strain,  and  fi{a,z,t)  is  the  internal 
shear  stress  acting  on  the  inner  cylinder  in  the  longitudinal  direction.  At  the  interface  of  the 
cylinders,  the  displacements  and  stresses  of  the  inner  cylinder  are  matched  to  the  displacements 
and  stresses  of  the  outer  cylinder.  These  boundary  interface  equations  are  written  as 


Ur(b_,z,t)  =  Ur(b+,z,t)  , 

(48) 

Uz(b_,z,t)  =  Uzib+,z,t)  , 

(49) 

Gpr{b—^Z,t)  Gp.(b^,Z,t)  , 

(50) 

OjrzQ^—iZyi)  —  Gf2{b^,z,t)  , 

(51) 

where  the  minus  subscript  corresponds  to  equations  evaluated  using  the  inner  cylinder 
properties  and  the  plus  subscript  corresponds  to  equations  evaluated  using  the  outer  cylinder 
properties.  The  normal  stress,  strain,  and  radial  forces  acting  on  the  outer  cylinder  at  radius  c 
are  related  by 

Crr{c,z,t)  =  (A2  +  2p.2)erric,2,t)  +  hieQQ{c,z,t)  +  Xiezz{c,z,t)  =  Poic,z,t)  ,  (52) 
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where  pQ(c,z,t)  is  the  external  pressure  acting  on  the  cylinder  in  the  radial  direction,  and  the 
subscript  2  denotes  material  properties  of  the  outer  cylinder.  The  shear  stress,  strain,  and 
longitudinal  forces  acting  on  the  outer  cylinder  at  radius  c  are  related  by 

(Jrz  (c,  z,  t)  =  2p2£rz  (c,  z,  0  =  fo  (c,  z,  t)  ,  (53) 


where  CTrzic,z,t)  is  the  shear  stress,  £rz{c,z,t)  is  the  shear  strain,  and  fo(c,z,t)  is  the  external 
shear  stress  acting  on  the  outer  cylinder  in  the  longitudinal  direction. 

Combining  equations  (36),  (37),  (38),  (40),  (42),  and  (46)  yields  the  normal  stress  in 
terms  of  the  potential  functions  g  and  hg  atr  =  a  as 

dr^  a  dr  dr^ 

where  is  the  magnitude  of  the  normal  force  acting  on  the  interior  of  the  inner  cylinder. 
Combining  equations  (39),  (44),  (45),  and  (47)  yields  the  shear  stress  in  terms  of  the  potential 
functions  g  and  Ag  at  r  =  a.  The  result  is 

+  +  ,  (55) 

dr  \  a^J  a  dr  dr^ 

where  Fi  is  the  magnitude  of  the  externally  applied  shear  stress  acting  on  the  interior  of  the 

inner  cylinder.  Equations  (40)  and  (48)  are  combined  to  yield  the  radial  displacement  interface 

equation 

^i^-ikhe(b_)  =  ^^-ikh0(b+),  (56) 

and  equations  (41)  and  (49)  are  combined  to  produce  the  longitudinal  displacement  interface 
equation 

r  dr  r  dr 

Equations  (36),  (37),  (38),  (40),  (42),  and  (50)  are  combined  to  yield  the  normal  radial  stress 
interface  equation 
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(Ai  +  +  -  ^ik^g(b. )  -  l^iik  = 

dr  o  dr  dr 

(A;  +  2h2)^^^  +  -  hl^'^g(K)- 2H2lk  .  (58) 

dr  0  dr  dr 

and  equations  (39),  (44),  (45),  and  (51)  are  combined  to  produce  the  shear  stress  interface 
equation 


lniik 


ndg(b-) 


a  /ii 


+  -7^  +  ^ 


dheib.)  .  ..  d%(b.) 


b  dr 


Ifilik 


+  ib+ )  + 


/^2  dh^jb^)  .  ..  d%{b^) 


b  dr 


+  (59) 

dr^ 


Combining  equations  (36),  (37),  (38),  (40),  (42),  and  (52)  yields  the  normal  stress  in  terms  of 
the  potential  functions  g  and  at  r  =  c  as 

(A,  +  P.  .  (60) 

dr  c  dr  dr 

where  is  the  magnitude  of  the  normal  force  acting  on  the  exterior  of  the  outer  cylinder. 
Combining  equations  (39),  (44),  (45),  and  (53)  yields  the  shear  stress  in  terms  of  the  potential 
functions  g  and  at  r-c .  The  result  is 

'^rr  I  (n  .  /^2  dh^ic)  ,  „  d\(c)  _  ^ 

^^^2ik—  +  ]^^i2k—^jheic)  +  - — ^  +  ^^2—^  =  Fo^  (61) 


where  Fq  is  the  magnitude  of  the  externally  applied  shear  stress  acting  on  the  exterior  of  the 
outer  cylinder.  Implicit  in  equations  (54)-(61)  is  the  assumption  that  the  external  loads  on  the 
cylinder  are  occurring  at  a  definite  frequency  and  wavenumber. 

Inserting  equations  (30)  and  (33)  into  equations  (54)-(6I),  applying  the  recurrence 
relationships  of  the  first-  and  second-kind  standard  Bessel  functions,  and  then  rewriting  as  an 
eight-by-eight  system  of  linear  equations  results  in 
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0 

^3 
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where  the  matrix  coefficients  are  given  as 

ail  -  j/iCaja)  ,  (63) 

«12  =[(-Ai -2;ii)af  - Ai/:2]ro(aia)  +  (^^^yi(aia)  ,  (64) 

ai3  =  (-2/ii/^/3i)7o(A«)  +  [^]/i(Afl)  .  (65) 

ai4  =  (-2Mii/c/3i)yo(/^a)  +  ,  (66) 

fll5  =  0,  (67) 

ai6=0>  (68) 

ai7=0,  (69) 

ai8=0.  (70) 

021  =  (-2^li/cai)Ji(aia)  ,  (71) 

a22=(-2Mii^oci)ri(aia)  ,  (72) 

a23=[^l(^^-A^)]7i(i3i«)  .  (73) 

a24=lMia^-fil)Wl(M  ,  (74) 

^25=0,  (75) 
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o 

II 

(76) 

It 

o 

(77) 

II 

o 

(78) 

a2,\  =  -aiJi{aib)  , 

(79) 

032  =  -aiY\{aib) , 

(80) 

033  =  -ikJii^b)  , 

(81) 

fl34=-»^n(W  . 

(82) 

^35  =  CC2Ji{a2b)  , 

(83) 

^36  =  > 

(84) 

a2i=ikJ\{fi2b)  , 

(85) 

a-ii^ikY^ihb)  , 

(86) 

a/^i=ikJQ(aib)  , 

(87) 

042  =  ikY^iciib)  , 

(88) 

^^43  =  A-^oCA^) . 

(89) 

fl44  =  AJ^o(A^)  . 

(90) 

045  =  -ikJo(a2b)  , 

(91) 

«46  = -'*^o(«2^)  . 

(92) 

^47  =  ~A‘^o(A^)  • 

(93) 

fl48  =  -AJb(A6) . 

(94) 

O51  =[(-Ai  -2/li)af  _ 

\h^cc\b)  . 

(95) 

052  =  [m,  - 

|li(«i*)  . 

(96) 
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a53  =(-2/Zii/:i8i)/oW)  +  f^Vl(A^)  ’  (97) 

\  b  J 

<i54=(-2/iittA)l'o(A*)  +  f^yi(Af>)  ,  (98) 

055  =[(A2  +2/12)02  +^*:^]^(a2*)-(^^Vl(«2'’)  .  (99) 

a56  =  [(^  +  2M2)a|  +  22t^])'o(a2*)-[^^yi(“2«  .  (100) 

oj,  =  (2/l2ilrft,)yo(ft«-f^^^Vl<fcW  ■  (*0') 

\  b  J 

fl58  =  i2H2ikP2)Yo(P2l>)  -  .  (102) 

%i=(-2^iikai)Jiiaib)  ,  (103) 

062  =(-2/Xi//:ai)}i(aiZ>)  ,  (104) 

a63=[^ii(k^-pi)]Jlim  ^  (105) 

064  =  \-^llik^  -  P^)]YliPlb)  ,  (106) 

065  =(2/i2i/^a2)*^l(«2^)  .  (107) 

066=(2/i2i^a2)li(«2*)  .  (108) 

067  =  ,  (109) 

a^=l-li2ik^-Ph]Yiim  ,  (110) 

a7i=0.  (Ill) 

^72=0,  (112) 

a73=0,  (113) 

^74=0,  (114) 

^75  =[(-^  -2/i2)a2  -^2^^]-/o(a2^^)  +  (^^]7i(a2c)  ,  (115) 
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“If,  =  Y-h  - 2/i2)«2  - , 

(116) 

a77  =  (-2/t2//:)52)/o()32c)  +  (^^^j7l()32c)  . 

(117) 

fl78  =  (-2/i2//:i52)lb(i32C)  +  f^^Vi(i%c)  , 

\  c  J 

(118) 

o 

II 

(119) 

fl82  =  0  , 

(120) 

^83  =  0  » 

(121) 

084=0  . 

(122) 

a^^  =  (-2p2ikoc2)J\(oc2c)  , 

(123) 

^86  =  (“2/X2l^«2)ll(®2<^)  > 

(124) 

aiT  =  lH2(k^-^)Vi(P2c)  . 

(125) 

«88=[#'2«:^-j32)»'l(ftc)  . 

(126) 

The  interior  pressure  field  at  the  inner  cylinder  a  is  now  broken  into  two  components: 
one  models  the  normal  wall  pressure  in  the  turbulent  boundary  layer  and  the  other  models  the 
acoustic  pressure  in  the  fluid.  This  interior  pressure  term  is  written  as 


Piia,z,t)  =  Pihi(a,z,t)- p(a,z,t)  ,  (127) 

where  the  subscript  tbl  denotes  the  term  associated  with  the  turbulent  boundary  layer  and  is 
discussed  in  the  next  section.  The  acoustic  pressure  in  the  fluid  is  derived  from  the  wave 
equation  written  in  cylindrical  coordinates  as 


p(r,z,t)  _  2 


d^p(r,z,t)  ^  1  dp(r,z,t)  ^  d^p(r,z,t) 


dr^ 


dr 


dz^ 


(128) 
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where  Cy  is  the  compressional  wave  speed  in  the  fluid  (m/s).  The  assumption  is  made  in 
equation  (128)  that  Mach  number  effects  on  the  wave  equation  are  negligible  because  the  fluid 
is  flowing  at  a  sufficiently  slow  speed.  The  boundary  condition  on  the  fluid/structure  interface 
is  written,  using  conservation  of  linear  momentum,  as 


dpia,z,t)  _  ^  d\ia,z,t) 
dr  -  *2 


(129) 


The  acoustic  pressure  field  is  modeled  as  a  magnitude  with  definite  frequency  and  wavenumber 
content  as 


p{r,z,t)  =  P{r)e^^e^^ 


Inserting  equation  (130)  into  equation  (128)  yields 


d^P{r)  ^  1  dP{r)  ^ 
dr^  r  dr 

The  solution  to  equation  (131)  is 


r 


Cf 


P(r)  =  0  . 


(130) 


(131) 


P{r)  =  C^jQ{yr)  , 

where  Cg  is  a  constant  (shown  implicitly  in  equation  (135))  and 


(132) 


y  = 


(133) 


It  is  implicit  in  the  solution  to  equation  (131)  that  the  pressure  in  the  interior  fluid  is  finite. 
Writing  the  radial  displacement  as 


Uj.{r,z,t)  =  Uj.{r)e^^e^^ 


(134) 


and  then  inserting  equation  (132)  into  the  boundary  condition  (equation  (129))  yields 

=  ,  (.35) 

r/i()o) 

It  is  now  noted,  based  on  equations  (12),  (30),  (33),  (40),  and  (134),  that 

Uria)  =  -CiaiJiiaia)  - CiaiY^ia^a) - - C^ikYi^^a)  .  (136) 
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Combining  equations  (135)  and  (136)  yields 


o/  X  ^  P/fi)Vo()c)ai/i(aia)  pfOp-jQ{ya)aiYi{aia) 

P(a)  =  Cl  — - +  Co  — - 

rM'io)  yfiiio) 

PfO^Joi'iaWii^a)  pfQ)^Jo('p)ikYi(^a) 

+  Co  — - +  Ca  — - 

x^i()o)  yf\{y3) 


(137) 


Equation  (137)  is  inserted  into  equation  (62),  which  results  in  the  A  matrix  coefficients  a^, 
a^2>  ^i3»  ^14 modified  to 


an  =  [(--^1  -  2pi)af  -  Xik^^Joiaia)  ■ 
a\2  -  [(““^1  “  2/Xi)a^  -  jyo(«ifl)  • 


a 


yfiiyi) 


2pi«i  ^  PfGrJQ(ya)a\ 
yJiiya) 


a 


Y\ia\a) 


(138) 


(139) 


ai3  =  i-2piik^)J Q{^a)  + 


^Ipiik  ^  PfaP-JQ{ya)ik^ 


a 


yJ\{ya) 


hiM . 


J 


and 


(140) 


^14  =  (•"2pij/:^)}o(A^)  + 


^Ipiik  ^  pfOp'JQ{ya)ik^ 


V 


a 


YiiM  . 


(141) 


^i()o) 

✓ 

The  remaining  terms  in  the  A  matrix  are  unchanged.  Additionally,  the  external  shear  stresses  at 
a  and  c  and  the  external  normal  pressure  at  c  are  assumed  to  be  equal  to  zero,  resulting  in  a  load 
vector  of 


Ptbl 

0 

0 

0 

0 

0 

0 

0 


(142) 
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Now  that  the  A  matrix  and  the  load  vector  are  known,  a  solution  to  the  constants  Cj  through 
C3  can  be  determined  by  inverting  the  A  matrix  and  multiplying  it  by  the  load  vector  in 
equation  (142).  Once  the  constants  are  known,  the  radial  displacement  at  c  divided  by  the 
interior  turbulent  boundary  layer  pressure  can  be  determined  with  an  expression  similar  to 
equation  (136): 

nk^coo)  =  =  -C5a2/i(a2^^)  -  Qa2>l(«2c)  - C-iikJii^c)  -  C^ikYii^c)  .  (143) 

nbi 

Figtire  7  is  a  plot  of  the  transfer  function  of  radial  displacement  divided  by  the  turbulent 
boundary  layer  pressure  versus  wavenumber.  The  frequency  of  the  transfer  function  is  250 
Hz,  the  inner  radius  (a)  is  0.00159  m  (1/16  in.),  the  middle  radius  {b)  is  0.00238  m  (3/32  in.), 
and  the  outer  radius  (c)  is  0.00191  m  (3/4  in.).  For  the  inner  cylinder,  the  dilatational 
(compressional)  modulus  (£j)  is  (l+0.18i)  x  10^  N/m^,  the  Poisson's  ratio  (u^)  is  0.49 
(dimensionless),  the  shear  modulus  (Gj)  is  3.4(l+0.18i)  x  10^  N/m^,  and  the  density  is  1000 
kg/m^.  For  the  outer  cylinder,  the  dilatational  modulus  (£2)  is  1.9(1+0.301)  x  10^  N/m^,  the 
Poisson's  ratio  (1)2)  is  0.49  (dimensionless),  the  shear  modulus  is  6.4(1+0.301)  x  10^ 
N/m^,  and  the  density  is  930  kg/m^.  The  speed  of  a  compressional  wave  in  the  fluid  (Cy)  is 
1500  m/s  and  the  density  of  the  fluid  (py)  is  1000  kg/m^.  Note  that  the  dynamic  range  of  the 
magnitude  of  the  transfer  function  is  approximately  100  dB  in  the  wavenumber  range  -750  <  k 
<  750  rad/m.  The  relative  minima  in  the  transfer  function  that  occurs  at  ±64  rad/m  correspond 
to  a  fluid/structure  interaction  bulge  wave.  The  maxima  of  the  transfer  function  correspond  to 
shear  wave  propagation. 
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2.4.  THE  SOURCE 


The  source  of  energy  for  the  excitation  of  the  viscoelastic  material  under  consideration 
is  from  the  flow  of  fluid  through  the  small  diameter  flexible  cylinder.  It  is  assumed  that  the 
flow  through  the  pipe  is  fully  developed,  turbulent,  and  incompressible.  The  pressure  at  the 
fluid-solid  interface  is  known  as  the  “wall”  pressure.  The  space-time  characteristics  of  this 
energy  are  used  as  input  to  the  stmctural  response  model  for  the  viscoelastic  material. 

No  measurements  of  the  wall  pressure  for  turbulent  flow  inside  a  small  diameter 
cylinder  are  known  to  exist.  Spectral  measurements  in  larger  cylinders,  however,  have  been 
found  to  agree  well  with  those  made  on  a  flat  plate  (Bakewell  et  al.,  1962).  The  source  model 
of  the  normal  pressure  of  the  turbulent  boundary  layer  in  the  cylinder  is  therefore  derived  from 
a  previously  formulated  flat  plate  wall  pressure  spectrum  model  (Chase,  1987).  Because  it  is 
convenient  to  work  in  the  spectral,  or  Fourier-transform,  domain,  the  model  represents  the 
two-wavenumber/frequency  spectrum  of  the  wall  pressure  as  opposed  to  the  space-time 
correlation  function.  The  spanwise  dimension  of  the  flat  plate  is  replaced  by  the  circumferential 
dimension  of  the  cylinder  in  which  the  fluid  flows.  This  translation  is  achieved  by  applying  a 
filter  function  in  the  spanwise  wavenumber  domain  to  the  two-wavenumber/frequency 
spectrum.  The  result  is  the  circumferentially  averaged  or  “mode  0”  component  of  the 
streamwise  wavenumber-frequency  spectrum  of  the  pressure  inside  the  cylinder. 

The  model  of  the  two-wavenumber/frequency  spectrum  of  wall  pressure  for 
incompressible  flow  over  a  flat  plate  is  given  as  (Chase,  1987) 


\P(k.k2.COof  -  f  ,512  \cT[k^*li) 


[K. 


l  +  {b5) 


Kl-¥{b8) 


-2 


+  ^2  +  ibSy 


+  Cif4k 


(144) 


where  k  is  the  wavenumber  of  the  streamwise  component  of  the  flow  (rad/m);  ^2  is  the 
wavenumber  of  the  spanwise  component  of  the  flow  (rad/m);  v*  is  the  friction  velocity  (m/s); 
5  is  the  boundary  layer  thickness  (m);  Cj-,  C^,  and  b  are  constants;  and 


^2 

(hv*r 


+  k2  , 


(145) 
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where  is  the  convection  velocity  (m/s)  and  /i  is  a  constant.  A  positive  value  of  k 
corresponds  to  fluid  energy  moving  downstream.  The  convection  velocity  is  assumed  to  be 
equal  to  the  mean  flow  velocity,  u,  and  they  will  be  used  interchangeably  in  the  following 
derivations.  The  friction  velocity  by  definition  is 


I 

I 


(146) 


where  is  the  shear  stress  at  the  wall  (N/m^).  The  mean  wall  shear  stress  inside  a  cylinder 
is  approximated  by 


^  -Pf^ 

''W  ~  o  ’ 


(147) 


where  /  is  the  friction  factor.  For  smooth  pipes,  the  friction  factor  is  (Fox  and  McDonald, 
1985) 


/  = 


0.3164 


(Red) 


0.25 


(148) 


where  Re^j  is  the  Reynolds  number  of  the  flow  calculated  using  the  maximum  (centerline)  pipe 
velocity.  The  ratio  of  the  centerline  velocity  to  the  mean  flow  is  (Fox  and  McDonald,  1985) 


IrP' 

U~  (n  +  l)(2n  +  l)  ’ 


(149) 


where  U  is  the  centerline  velocity.  The  parameter  n  is  the  inverse  of  the  exponent  used  in  the 
power  law  velocity  profile,  which  is  often  used  to  model  velocity  profiles  in  turbulent  flow. 
For  the  range  of  Reynolds  numbers  in  this  study,  a  value  of  n  =  4  is  applicable. 

For  fully  developed  pipe  flow,  the  boundary  layer  thickness  S  is  equal  to  the  cylinder 
inner  radius.  To  effectively  transform  the  flow  over  a  flat  plate  into  flow  inside  a  cylinder,  a 
filter  function  is  applied  to  the  two-wavenumber/frequency  spectrum  (equation  (144)).  The 
filter  function  is  that  of  a  simple  linear  continuous  sensor  with  a  length  equal  to  the  inner 
circumference  (2kS)  of  the  cylinder.  The  integral  in  spanwise  wavenumber  of  the  product  of 
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the  filter  function  and  the  two-wavenumber/frequency  spectrum  provides  a  new  model  of  the 
mode  0  component  of  the  streamwise  wavenumber-frequency  spectrum  and  is  given  as 


oo 

=  J  \P{k,k2,(OQf 


sin(^2^^) 

kind 


1 

dk2  . 


(150) 


Figure  8  is  a  plot  of  the  wavenumber-frequency  spectrum  of  mode  0  in  a  pipe.  The  top 
figure  is  plotted  using  a  wavenumber  range  of  ±750  rad/m,  and  the  bottom  figure  is  plotted 
using  an  expanded  wavenumber  range  of  ±5(X)0  rad/m  to  show  the  effects  of  convective  ridge 
energy.  The  radius  of  the  pipe  is  0.00159  m  (1/16  in.),  the  mean  flow  is  0.535  m/s,  the 
density  of  the  fluid  is  1000  kg/m^,  the  kinematic  viscosity  is  1.01  x  10*^  ra^/s,  and  the 
frequency  is  250  Hz.  The  constant  his3,  b  is  0.75,  Cj  is  0.00467,  and  Cj^  is  0.155.  The 
transfer  function  (equation  (150))  in  figure  8  is  multiplied  by  2;r  to  convert  the  frequency 
portion  of  the  units  from  l/(rad/s)  to  1/Hz.  Although  the  energy  level  of  the  function  reaches  a 
maximum  at  the  convective  wavenumber,  this  energy  is  filtered  out  of  the  system  because  of 
the  transfer  function  of  the  combined  fluid  and  solid  structure.  The  dominant  measured  effects 
of  the  fluid  occur  between  -250  and  250  rad/m,  which  is  the  unaliased  wavenumber  region  of 
the  array. 


29 


TR  11,073 


3.  EXPERIMENT 

The  experiment  was  conducted  in  an  anechoic  chamber  with  the  laboratory 
configuration  and  array  previously  described.  The  data  were  collected  in  the  time  domain  using 
a  Labview  Data  Acquisition  system  residing  on  a  Macintosh  computer.  The  time  domain  data 
were  processed  as  96  ensemble  averages,  consisting  of  1024-point  nonoverlapping  fast  Fourier 
transforms,  with  a  Hanning  window  used  for  conversion  to  the  frequency  domain.  The  spatial 
domain  data  were  processed  as  a  32-point  fast  Fourier  transform,  with  three  zero  pads  to  every 
one  data  point.  This  approach  provides  a  convenient  method  of  interpolating  the  spatial  data  in 
wavenumber. 

The  surgical  tube  inner  radius  (a)  is  0.00159  m  (1/16  in.),  and  outer  radius  (b)  is 
0.00238  m  (3/32  in.).  The  distance  from  the  middle  of  the  surgical  tube  to  the  surface  of  the 
urethane  (c)  is  0.00191  m  (3/4  in.).  The  material  properties  of  the  surgical  tube  and  block  were 
measured  on  a  Metravib  viscoanalyseur  and  are  listed  in  tables  1  and  2,  respectively.  The 
water  had  a  density  of  1000  kg/m^,  a  compressional  wave  speed  of  1500  nVs,  and  a  kinematic 
viscosity  of  1.01  x  10"^  m^/s.  The  four  terms  of  equation  (1)  were  assembled  and  numerically 
integrated  from  -750  to  750  rad/m.  Although  the  limits  of  integration  in  equation  (1)  are  from 
minus  infinity  to  plus  infinity,  the  wavenumber  filter  rejects  most  of  the  high  wavenumber 
energy  (e.g.,  \k\>Q)/U(.)  and  thus  allows  only  low  wavenumber  energy  to  be  passed  to  the 
sensors. 

Because  the  measured  displacements  are  extremely  small,  a  reliable  method  to  calibrate 
the  array  has  not  yet  been  developed.  When  the  theoretical  results  are  compared  to  the 
experimental  results,  the  maximum  value  (in  wavenumber)  of  the  theoretical  results  is  matched 
to  the  maximum  value  of  the  experimental  results.  Figures  9  to  15  show  these  results  plotted  at 
100,  150,  200,  250,  300,  350,  and  400  Hz,  respectively.  The  theory  is  plotted  with  a  solid 
line  and  the  experimental  data  with  X’s.  The  mean  flow  in  the  experiment  was  0.535  m/s,  with 
a  corresponding  diameter-based  Reynolds  number  of  1680.  Although  this  is  slightly  below  a 
fully  developed  turbulent  flow  Reynolds  number  of  2000,  it  is  close  to  the  higher  Reynolds 
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number  of  1500  that  has  been  observed  downstream  of  stenosis  in  laboratory  experiments 
(Jones  and  Fronek,  1988). 


Table  1.  Material  Properties  of  the  Surgical  Tube 


Compressional  Modulus 

Shear  Modulus 

Frequency  (Hz) 

Real  (£)  N/m^ 

Imag  (E)  N/np- 

Real  (G)  N/trP 

Imag  (G)  N/m^ 

100 

10.  X  10^ 

1.5  X  10^ 

3.4  X  10^ 

0.50x  105 

150 

10.  X  105 

1.6  X  10^ 

3.4  X  105 

0.54  X  10^ 

200 

O 

X 

d 

1.7  X  10^ 

3.4  X  10^ 

0.57  X  105 

250 

10.  X  10^ 

1.8  X  10^ 

3.4  X  10^ 

0.60  X  10^ 

300 

10.  X  lO^ 

1.9  X  10^ 

3.4  X  lO^ 

0.64  X  10^ 

350 

10.  X  10^ 

2.0  X  10^ 

3.4  X  10^ 

0.67  X  10^ 

400 

p 

X 

o 

2.0  X  10^ 

3.4  X  10^ 

0.67  X  10^ 

Table  2.  Material  Properties  of  the  Urethane 


Compressional  Modulus 

Shear  Modulus 

Frequency  (Hz) 

Real  (E)  N/m^ 

Imag  (E)  N/nP 

Real  (G)  WnP 

Imag  (G)  N/m^ 

100 

2.1  X  10^ 

0.36  X  105 

0.70  X  105 

0.12  X  105 

150 

2.1  X  10^ 

0.48  X  105 

0.70  X  105 

0.16  X  10^ 

200 

2.1  X  10^ 

0.59  X  10^ 

0.70  X  105 

0.20  X  105 

250 

1.9  X  10^ 

0.57  X  10^ 

0.64  xl05 

0.19  X  105 

300 

1.9x  105 

0.67  X  105 

0.64  X  105 

0.22  X  105 

350 

1.8  X  10^ 

0.68  X  105 

0.60  X  105 

0.23  X  105 

400 

1.7  X  105 

0.71  X  105 

0.57  X  105 

0.24  X  105 
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4.  DISCUSSION 

The  turbulent  energy  present  in  the  embedded  tube  was  measured  at  the  surface  of  the 
solid  with  the  array.  The  energy  traveled  to  the  surface  by  means  of  a  propagating  shear  wave 
that  originated  at  the  tube  wall.  The  effects  of  compressional  wave  motion  were  small,  and 
these  waves  did  not  show  up  as  maxima  in  the  wavenumber  cuts.  The  model  was  rerun  with 
varying  compressional  wave  speeds  that  produced  results  similar  to  the  original  results. 
Changing  the  stiffness  of  the  embedded  tube  had  a  moderate  effect  on  the  dynamics  of  the 
modeled  system.  Changing  the  modeled  shear  modulus  of  the  embedded  tube  and/or  the  solid 
moved  the  location  of  the  peaks  on  the  wavenumber  plots.  This  factor  must  be  considered 
when  the  array  is  used  to  make  in  vitro  measurements  on  human  arteries. 
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The  shapes  of  all  the  theoretical  frequency  plots  are  in  relative  agreement  with  the 
experimental  results.  The  model  accurately  predicts  the  shear  wave  characteristics  that  are 
present  in  the  solid.  In  figures  10,  11,  and  12,  the  magnitude  of  the  near-zero  wavenumber 
regions  do  not  agree  exacdy  in  shape,  which  is  probably  due  to  the  curved  surface  assumption 
used  to  approximate  the  flat  surface  on  the  urethane  block.  At  low  wavenumbers,  the 
wavelengths  are  long  and  the  effect  of  surface  tension  is  likely  to  be  more  significant  than  at 
higher  wavenumbers.  This  magnitude  mismatch  could  also  be  the  result  of  sensor-to-sensor 
magnitude  and  phase  mismatch.  Additionally,  there  is  some  high  wavenumber  mismatch 
between  the  data  and  the  model,  which  could  be  due  to  the  phase  centers  of  the  elements  of  the 
array  not  being  exactly  0.5  inches  apart.  Slight  phase  center  mismatch  will  result  in 
experimental  results  that  are  shifted  slightly  in  wavenumber. 

The  source  model  assumed  fully  developed  turbulent  flow  in  a  long  pipe  and  is  derived 
from  a  fiat  plate  model  of  the  wall  pressure  spectrum.  The  accuracy  of  results  based  on  these 
premises  requires  quantification.  Effects  of  transitional  (laminar  to  turbulent)  flow, 
nonhomogeneous  tubes  (e.g.,  one  with  a  stenosis),  and  small  radius  (curvature)  effects  should 
also  be  examined  more  closely  in  further  expaimentation. 

This  research  illustrates  the  need  to  accurately  calibrate  the  sensors.  Because  the  peak 
value  of  the  measured  response  was  matched  to  the  peak  value  of  the  theoretical  response,  the 
comparison  of  the  experimental  and  theoretical  results  is  relative.  Although  this  method 
provides  a  good  agreement  across  wavenumber,  a  better  approach  would  be  to  make  a 
calibrated  measiuement  and  then  compare  the  measured  response  to  the  theoretical  response. 

The  wavenumber-frequency  beamformer  used  in  this  research  was  extremely  accurate 
in  resolving  the  wavenumber  content  of  the  space-time  field.  However,  it  is  likely  that  a  high 
resolution  (adaptive  or  optimal)  beamformer  would  provide  a  better  signal  processing  method 
to  locate  spatial  sources.  Such  a  beamformer  would  incorporate  the  spherical  characteristics  of 
the  wave  and  a  corresponding  loss  term. 
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5.  CONCLUSIONS 


The  space-time  field  of  flow  through  a  long  pipe  surrounded  by  an  elastic  medium  was 
modeled  and  experimentally  verified.  A  comparison  of  the  model  to  the  experiment  showed 
broad-based  agreement  in  wavenumber  and  frequency. 
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